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Abstract 


In response to a systematic methodology assess- 
ment program directed to the aeroelastic stability 
of hingeless helicopter rotor blades, improved basic 
aeroelastic reformulations and new formulations 
relating to structural sweep have been achieved. 
Correlational results are presented showing the 
substantially improved performance of the G400 aero- 
elastic analysis incorporating these new formula- 
tions. The formulations pertain partly to sundry 
new solutions to classic problem areas, relating to 
dynamic inflow with vortex-ring state operation and 
basic blade kinematics, but mostly to improved 
physical modeling of elastic axis offset (structural 
sweep) in the presence of nonlinear structural twist. 
Specific issues addressed are an alternate modeling 
of the AEI torsional excitation due to compound 
bending using a force integration approach, and the 
detailed kinematic representation of an elastically 
deflected point mass of a beam with both structural 
sweep and nonlinear twist. 
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Component of load distribution in radial (x ) 
direction, ND 

Tension at blade section, or rotor thrust, 
as appropriate, lbf . 

Coordinate transformation matrix relating "5" 
and "6" coordinate systems, due to structural 
sweep , ND 

Inward radial (x^) foreshortening of blade 
element point due to combination of built-in 
sweep and elastic deformation, ND 
inflow parameter 

Elastic deflections in the edgewise and flatwise 
directions, respectively, ND 

Uniform component of momentum induced velocity, ND 

Cosine and sine components, respectively, of 
momentum induced velocity, ND 

Deflection correction functions due to first order 
twist effects, ND 

Deflection correction terms due to second order 
twist effects, ND 

Components of position vector in the "5" system 
(rotating, coned and lagged), ND 
Built-in offset distances of elastic axis from 
*5 axis in inplane and out-of-plane directions, 
respectively, ND 

Built-in offset distances of elastic axis from 
x 5 axis, in edgewise and flatwise directions, 
respectively, ND 
Built-in blade precone, deg. 
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Tip loss factor 

Rotor roll and pitch moment coefficients, 
respectively, (moment / c ^r.“R 5 ) , ND 
Rotor thrust coefficient (T/pr^R**), ND 

Section bending stiffness in flatwise and edgewise 
directions, respectively, lb-in' or ND, 

(AEI-EI -El ) 
z y 

Induced velocity function, ND 
Tension cosine resolution function, ND 

Induced velocity gradient factor, ND 

Section shear load distributions in directions 

of "5" coordinate system, ND 

Blade k'th edgewise modal response variable 

Blade i'th flatwise modal response variable 

Section moment load distributions about axes 

in the ”5" coordinate system, ND 

Blade j ' th torsion modal response variable 

Rotor radius, ft. 

Blade spanwise coordinate, measured from offset 
in x^ direction, ND 
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Built-in precone outboard of pitch bearing 
(negative droop), deg. 

Nonlinear j'th torsion modal weighting functions 
for torsion excitation due to edgewise and flat- 
wise force loadings, respectively, ND 

Nonlinear j 1 th torsion modal weighting functions 
for torsion excitation due to flatwise and edge- 
wise moment loadings, respectively, ND 
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Inplane and out-of-plane slope projection angles, 
respectively, defining blade element orientation, 

rad . 

Deflection mode shape for the k’th edgewise 
normal mode, ND 

Deflection mode shape for the i'th flatwise 
normal mode, ND 

Deflection mode shape for the j'th torsion normal 
mode, ND 

Total local blade pitch angle, radians 
Elastic torsion deflection angle, radians 

Collective pitch angle, deg. 


A 



Structural sweep angle projection onto x -y 
plane, rad. 5 3 

Structural sweep angle projection onto x -z 
plane, rad. 5 5 
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Inflow ratio with spanwise and azimuthal 
variability, ND 

Normalized rotor through flow parameter, ND 
Part of uniform inflow arising from rotor 
forward flight, ND 

Uniform component of variable inflow, ND 

Rotor advance ratio, ND 
Air density, lb-sec 2 /ft^ 

Alternately, rotor solidity, and real part of 

eigenvalue, ND 

Blade azimuth angle, rad. 

Rotor rotation speed, rpm. 


Subscripts and Superscripts 


< >e 

( * ^ 
( * ) 

( )' 

( A ) 


Due to elastic deformation 
Defined at the elastic axis 
Differentiation with respect to ^ 
Differentiation with respect to (r/R) 

Denotes ^valuation at zero collective 
angle as applied to deflections 

Introduction 


For most production helicopter design applica- 
tions, the principal role of contemporary comprehen- 
sive rotor aeroelastic analyses has been that of 
providing calculations of forced structural responses 
and, in particular, of blade dynamic stresses. The 
United Technologies Corporation family of G400 rotor 
aeroelastic analyses comprises such a comprehensive 
analysis technology and has undergone extensive 
development in the last ten years with this principal 
role as a prime objective. The present G400 techno- 
logy has evolved from an analysis originally formu- 
lated for the unique aeroelastic characteristics of 
the composite bearingless rotor. That analysis 
represented an advancement in the state-of-the-art 
with regard to the modeling of rotors with time- 
variable, nonlinear structural twist and multiple 
structural redundancy, as described in Reference 1. 
The G400 technology which has evolved now includes a 
family of four actively used versions with a 
completely general range of applicability in rotor 
type (articulared, hingeless, teetered and gimballed) 
and vehicle application (helicopters, propellers and 
wind turbines). The mathematical modeling capabili- 
ties of the G400 analyses are summarized in Figure 1. 
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G400 AEROELASTIC ANALYSIS 


• BEAM BENDING AND TORSION MODES 

• STRUCTURAL SWEEP AND TWIST 

• UNSTEADY AIRLOADS 



COUPLED (LINEAR) TRANSIENTS, PERFORMANCE HARMONIC 

MODES AEROELASTIC (NONLINEAR) RESPONSES. 

AND STABILITY AEROELASTIC STABILITY, STRESSES 

FREQUENCIES CONTROL INPUTS 

Fig. 1 - Basic capabilities of G400 Aeroelastic 
Analyses . 


Of the two major solution types, eigensolution and 
time-history solution, the latter contains the most 
complete physical modeling of the blade aeroelas- 
ticity. This includes the dynamics, airloads, exci- 
tations and kinematic couplings with the full reten- 
tion of all nonlinearities which have been identified 
as being potentially germaine to the aeroelastics . 
Thus, no nonlinearities have been deleted from the 
time-history solution for reasons of mathematical 
convenience. Prior to 1983, the major documentation 
of the G400 technology was available only in Refer- 
ences 1 through 3. Since completion of the work 
reported herein, another major documentation source 
has become available (Reference 4). 

Within the context of only forced response 
calculations, limited harmonic response correlation 
studies have been performed. These have been 
conducted principally under corporate and contractual 
funding; References 5 and 6 are the available docu- 
mentations of this type of correlation study. 

Detailed aeroelastic stability correlation studies, 
however, had not been performed prior to the perfor- 
mance of the Integrated Technology Rotor/Flight 
Research Rotor (ITR/FRR) Methodology Assessment study 
(Reference 7) . One reason for the lack of G400 
stability correlation calculations is clearly the 
emphasis placed on forced response loads calculations 
by the principal users of the code. Another more 
logistical reason, however, is that over most of its 
development life the G400 analysis has been princi- 
pally a time-history solution analysis. As a result, 
the eigensolution capability had not kept pace with 
the increased sophistication of this time-history 
solution capability. Consequently, accurate stabi- 
lity calculations have typically required the use of 
transient time-history calculations. Such calcula- 
tions are generally both time and cost intensive and, 
hence, had been eschewed. Despite the cost disadvan- 
tage, however, time-history solutions present a 
distinct advantage in the calculation of transient 
stability, as is discussed in greater detail in a 
subsequent section. 

Under contract NAS2-10864, the in-house heli- 
copter version of G400 was exercised for stability 
correlation as part of this methodology assessment 
study. Initial results of this study were generally 
poor. The G400 stability predictions were deemed 
unacceptably inaccurate and a concerted corporate- 
sponsored methodology improvement project was 
initiated. The general results of this improvement 
project were completely successful. The stability 
predictive capability of G400 was definitely raised 
to an acceptably accurate level (giving good to 
excellent correlation results) while retaining a 
valid, mathematically consistent formulation. Over 
and above this immediate positive result, however, 
this methodology improvement study produced new 
formulations and revised existing ones; these 
formulations are of interest in their own right. 
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The nature of the detailed reformulations were 
of three main types: The first consisted of the 

detection and correction of outright errors in the 
programmed implementations of the existing derived 
equations. The second consisted of a sundry class 
of modifications wherein established aeroelastic 
methodology was extended from the generally accepted 
norm. And the third consisted of an improved repre- 
sentation of structural sweep. A discussion of the 
first type of reformulation is clearly inappropriate 
for publication and is omitted from further discus- 
sion. The second and third types of ref ormulation, 
however, constitute new knowledge and form the basis 
of this paper. The remainder of this paper is 
divided into three main sections: (1) a review of 

the pertinent G400/ITR correlation results, (2) a 
description of the sundry modifications arising from 
enhanced reformulations of existing theory, and 
(3) a description of the new formulations relating 
to structural sweep. 


Review of Pertinent ITR Correlation Results 

The ITR Methodology Assessment Study, as 
defined in Reference 7, concentrated on the aero- 
elastic stability characteristics of hingeless and/ 
or bearingless rotors both in hub-fixed and hub- 
flexible configurations. Particular emphasis was 
placed on the stability of the already lightly damped 
blade edgewise (inplane) mode as affected by coup- 
lings with the blade flatwise (out-of-plane) and 
torsion modes, and with the flexible hub degrees-of- 
freedom. In all cases, the pertinent mode, whose 
stability characteristics were to be calculated, was 
characterized by relatively low reduced frequencies 
along the blade and for most conditions by an absence 
of stall. Hence, the stability phenomena could be 
assumed to be reasonably well-governed by conven- 
tional quasi-static airloads. 

The original results from applying G400 to the 
experimental correlational data were generally poor 
for most of the configurations defined in the study. 
Of particular significance were the poor correlations 
achieved with the simplest configuration: that of an 

isolated hingeless model rotor with no twist or 
cyclic pitch (configuration IIA, as described in 
detail in Reference 8). Although the other configu- 
rations were equally, if not more, important to the 
ITR study as a whole, only this configuration will be 
addressed in this paper because it was the primary 
vehicle which led to the enhancements to be 
discussed herein. 
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FLEXURE (STIFF OR SOFT IN TORSION) 


Fig. 2 Correlation cases for ITR configuration 
IIA, isolated hingeless rotor. 


The configuration IIA rotor stability data con- 
sisted of 6 distinct cases involving simple parameter 
variations in precone, 8^, droop, (-)AB, and torsional 
flexure stiffness, as shown in Figure 2. A measure of 
the torsional stiffness of the two flexures is 
a ff° r ded by the first torsional mode amplitudes near 
the blade root. For the stiff and soft flexures, the 
calculated torsion modal amplitudes (at the 3% span- 
wise location) were, respectively, .00013 and 0.12754 
For each of these parameter variations, the damping 
constant, o, was obtained as a function of blade 
collective angle, 6 0 , as shown in Figures 3a thru 3f. 
These figures present the experimentally obtained 
values together with the initial (12/81) G400 calcu- 
lations and the updated (5/83) ones. The improved 
correlation of the updated G400 results is apparent 
and is generally representative of all the results 

obtained by including the three types of reformula- 
tions. These figures will be referred to in the 
subsequent sections to illustrate the impact of the 
various specific reformulations. 


(a) Case 1 - stiff flexure, 6g= A3= 0. 



0 O . DEG 


Fig. 3 - Comparison of experimental results with 
initial and revised G400 calculations- 
conf iguration IIA. 
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(b) Case 2 - soft flexure, Bg= AB= 0. 






Sundry Reformulations Relating to Existing Theory 
Air Mass Dynamics 

Examination of Figures 3a and 3b shows an 
experimentally observable falloff of rotor stability 
at sufficiently high values of collective angle. 

The most obvious inaccuracy of the initial G400 
calculations is the premature falling off of the 
damping constant with increasing collective angle. 

The physics governing this attenuation of damping 
is twofold: First, increases in collective will 

necessarily increase the blade loading and, thereby, 
the static out-of-plane blade bending. This increase 
in static "bending will significantly impact on the 
effective pitch-edge coupling which, in large measure, 
defines the pitch-flap-lag stability. Secondly, 
increases in collective will also increase the pene— 
trat ion of the blade section angles-of-at tack into 
the near stall, high drag rise coefficient regime of 
airfoil operation. As shown in Reference 8, this 
regime of rotor operation'is generally destabilizing. 
The basic parameter common to and controlling each of 
these effects is the local blade section angle-of- 
attack. The angle-of-attack, however, is determined 
from both geometric and inflow contributions. From 
inspection of the initial G400 results it appeared 
that the section angle-of-attack vs. pitch angle 
relationship might be incorrect and such in fact 
was the case. 


and where v x , vj c and v lg are the uniform (zeroth 
harmonic), and first harmonic components of induced 
velocity, respectively. These components of induced 
velocity are governed by appropriate first order 
differential equations: 
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(4) 


where is a newly-defined rotor induced velocity 
function whose independent variable is taken to be 
the normalized through-flow parameter defined 
as follows: 


S, = sgn(x 0 C T )V^ + X z 0 / -J |C t |/2B 2 


(5) 


and the usual inflow parameter, v, is defined as: 


V 


M 2 + X 0 (Xo— V 0 ) 


(6) 


The G400 technology incorporates a representa- 
tion of air mass dynamics which closely conforms to 
the established state-of-the art (e.g. Reference 9). 
The major departure of the G400 technology from that 
typified by Reference 9 is twofold. First, the 
technology employs a nonperturbational, totally 
nonlinear form of the momentum equations. Second, 
in order to accommodate the high thrust loadings 
at which a wind turbine is capable of operating, the 
G400 technology employs an empirical correction 
procedure for simulating operation in or near the 
vortex-ring state. These ideas are summarized in 
the following development. The total (nonperturba- 
tional) form of momentum variable inflow is assumed 
to be as follows: 

X(r,^) = X ram - v 0 - r[(V| C + K v 0 ) cos ^ 4* v, $ sin^J (] 

where the Glauert factor, K, is approximated^ 10 ^ 
by the following simple expression: 

Cm /x 0 > 

Kr 3 I 2 +(m/X 0 ) (2i 

, where: 


X 0 z X RAM ~ V 0 (2b) 


For rotor operation well removed from the 
vortex ring state (|x|sl.4) the rotor induced 
velocity function, consists of two branches 
and is directly obtainable from standard momentum 
theory as given simply by 1/ | X | . For values of 
IM l e ss than 1,4 and especially approaching zero, 
the momentum representation breaks down (and 
eventually goes singular). Alternate empirical 
correction curves which connect the two valid 
momentum branches for values of X between -1.4 and 
+1.4 are suggested by material presented both by 
Gessow and Myers OD and by Lissaman, as shown 
in Figure 4. 



INFLOW PARAMETER, Sgn(XC T ) yx2^J /v /|^j7ip 


Fig. 4 Rotor induced velocity function. 
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Essentially the reformulations pertainxng to 
air mass dynamics which were included in the updated 
G400 technology were to include the signum function 
factor in the definition for X, as given in Equation 
(5) (in order to accommodate negative values of 
inflow, X), and to abandon the Lissaman data in favor 
of the Gessow and Myers data. For the configuration 
HA correlation cases, these changes resulted in 
values of :T which were above the momentum values 
compared with ones which were initially below, at 
the high thrust (high collective angle) conditions. 
Tliis correction to the formulation of the induced 
velocity function accounted for the difference in 
steady section angles-of-attack needed to bring the 
high collective pitch angle results into agreement 
with experiment. 

Basic Considerations of Blade Kinematics 

The high relative torsional stiffness of the 
"stiff" flexure, cases 3 and 5 of configuration IIA 
(see Figure 2), result in these cases taking on 
especially useful significance. For these two cases, 
the rotor blade is essentially rigid in torsion up to 
the point just outboard of the flexure. Thus, they 
are aeroelastically equivalent and should have the 
same stability characteristics. The experimental 
results shown in Figures 3c and 3e do confirm this 
supposition. 

Within the context of the G400 technology, 
however, cases IIA-3 and IIA-5 must be respectively 
modeled as a blade with a straight elastic axis 
preconed at a 5 degree angle, and as a blade without 
precone, but with a 5 degree bend in the elastic axis. 
The effective equivalency of cases IIA-3 and IIA-5 
thus forms the basis for validating the consistency 
of formulations especially with regard to elastic 
axis offset (structural sweep). 

The aeroelastic significance of both radial 
foreshortening and spanwise tension (treated in the 
subsequent subsections) is that they are each an 
important source of coupling between flatwise bending 
and edgewise bending. Because of the contributions 
of flatwise bending to radial foreshortening, flat- 
wise rate terms appear in the Coriolis force depen- 
dent terms in the edgewise equation. Similarly , 
because of the contribution of edgewise rate to the 
centrifugal force, edgewise rate terms appear in the 
tension terms in the flatwise equation. Because of 
the evident significance of these terms, a useful 
test for assessing the accuracy and self-consistency 
of the improved formulations was that the stability 
predictions for cases IIA-3 and IIA-5 be the same. 

Kinematics of Radial Foreshortening 

The original G400 development W invoked 
various principal assumptions which were intended 
to allow for advancement of the art of modeling 
nonlinear structural twist while avoiding unnecessary 
obfuscation caused by the inclusion of numerous 
nonlinear terms. Accordingly, the radial foreshor- 
tening of a mass element due to elastic bending, 


u g , was kept simplistic and assumed to be limited to 
that accruing from flatwise bending only. It was 
accordingly represented by a quadratic function in 
flatwise bending: 

u e= Tf m [/ 0 r <^i]qw,qw m 5 

(4) 

In the reformulated G400 technology , this 
restirctive assumption was relaxed. The two basic 
assumptions which were retained, expanded upon and- 
utilized as an alternative basis are as follows: 


1) The elastic (torsion axis is defined to be the 
spanwise locus of shear centers of the two- 
dimensional blade (beam) sections taken perpen- 
dicular to this spanwise locus. Note that this 
definition treats the elastic axis as an 
abstracted section property, as contrasted with 
what one would measure in a bench test of an 
actual curved beam. The built-in structural 
sweep (elastic axis offset), together with the 
elastic bending deflections, define an elastic 
axis which is generally a space-curve about which 
the local torsion deflection must take place. 

2) The arc length of the so-defined elastic axis is 
invariant both in toto and per blade segment. 
Radial foreshortening accrue entirely from the 
kinematics of bending and distributed torsion 
along the space-curve elastic axis. 

3) Local radial foreshortening is defined relative 
to the total extended arc length of the elastic 
axis. A hypothetical beam formed by the straigh- 
tening out of the arc length of the elastic axis 
and the elimination of all pitch and twist is 
herein defined to be the "equivalent beam." 


Contributions to radial foreshortening then 
accrue from (a) the built-in structural sweep, i.e. 
that which restores the equivalent beam to the origi- 
nal swept planform (b) first order (linear) functions 
of bending, arising from built-in structural sweep, 
(c) second order (nonlinear) functions of bending 
each with elastic torsion arising from built-in 
structural sweep, and (d) second order functions 
each of both flatwise and edgewise bending. 

These contributions are pictorially indicated in 


Figure 5. 
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Contributions to incremental raaiai tor 
shortening due to structural sweep and 
elastic deformations. 
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Each of these contributions can be 
modeled in a straightforward manner, and in lieu 
of the detailed development given in Reference 4, 
are simply stated as follows: 

(dAx), = dr-dx ' dr- ^/dr 2 - Ay| £A ~ Az| ea (3) 

* 1 — ^Ar)* — (A Zioea / Ar) 2 j 

r a / * / 

(dAx) 2 = dr[smAf 5 cosAe 5 z 6 e “ sinAe 5 cosAf 5 y6 c j (9) 

(dAx) 3 = cosA e 5o cosAf 5Q [i - ^/i - v e ' 2 - w* 2 ] dr 

= COS Ae 5o cosAf^ -^ ( v^ 2 + we' 2 )dr (10) 

where Ay^o^ and Azio^ are, respectively, the 
built-in changes per segment length of the chordwise 
and flatwise distances of the elastic axis from the 
reference, X5 , axis. And where A e ^ and Af^ are, 
respectively, the structural sweep angle projections 
onto the X 5 ~y 5 and X5-Z5 reference planes. 

The total elastic radial foreshortening at the 
center of the nth segment is then determined by the 
following integral: 

u en = [(dAx), + (dAx) 2 + (d Ax) 3 j (11) 


The details of this integration are straightforward 
but sufficiently tedious to be beyond the intent and 
usefulness of this paper. Symbolically, u e is 
finally given by: 

Ue = (DUE AO) + (DUEAFjJq*. + (DUEAE k ) • q VR 

+ (UELSET k ) q q 0 . + ( UELSFTj ) q* q 9 
J A J J 1 J 

+ i(UELASE km )q Vk q Vm + i(UELASF in )q W| q Wn (12) 

This formulation thus contains Equation 7 as a 
contributing term. 

Spanwise Tension Distribution 

Of all the terms appearing in the blade dynamic 
equations, the tension force is by far the greatest 
in magnitude and, by definition, qualifies as a 
"zeroth order” term. The difficulty in accurately 
modeling tension is that although it is a zeroth 
order term, the zeroth order component is equili- 
brated by other zeroth order effects (e.g. the steady 
blade airloads). Indeed, it can be well appreciated 
that the significant dynamics of rotor blades are 
determined by the higher order terms. Thus, even 


though tension is principally a zeroth order quantity, 
it still becomes important to model it with suffi- 
cient detail to capture the salient higher order 
effects . 

Tension has been typically calculated as the 
direct spanwise integration of the radial force 
loading, S x , outboard of the blade field point 
(center of blade segment). The radial force loading 
is, in turn, taken to be that due to centrifugal 
force and is thus dependent on the mass element 
radial position and inplane velocity, both of which 
include higher order terms. The formulations of 
the previous subsection, therefore, clearly impact 
on the calculation of centrifugal force. In addition 
to these reformulations, an additional higher order 
effect relating to tension was identified which 
subsequently led to the required self-consistency. 

In the reformulated G400 technology, account has 
been taken of the fact that tension is a vector 
whose local direction is determined by the orienta- 
tion of the beam element (blade segment). The cen- 
trifugal force on the other hand is a vector always 
oriented radially in the rotor rotation plane. 

Hence, tension and centrifugal force are not 
generally codirectional . Upon defining the out-of- 
plane and inplane projections of the skew angle, 
y, between these two vectors, as Yp an <3 Ye» 
respectively, the effect of non-codirectionality on 
tension can be written as: 

T(r) = F T (r) f f Sxjdr, (13) 

where : 

F T (r) = cosr(r)= ./-sin 2 r E -sin 2 r F (14) 

Reformulations Relating to Variable Elastic Axis 
Offset 

As originally formulated, the G400 technology 
assumed the elastic axis to define a space-curve as 
a result of combined flatwise and edgewise bending. 

In this case, the blade curvature is directly pro- 
portional to the elastic modal degrees-of-f reedom. 

Tli is situation consequently allowed for considerable 
simplicity in structural modeling especially with 
regard to the nonlinear torsion excitation resulting 
from combined flatwise and edgewise bending (the 
AEI term) For the case of built-in variable elastic 
offset (structural sweep) the accurate definition of 
such sweep in terms of its curvature becomes 
impractical. Also, while an approximation to the 
blade kinematics resulting from "small sweep" could 
be obtained heuristically by considering the struc- 
tural sweep to consist of "pre-bends" in the elastic 
axis, this procedure becomes suspect at moderate to 
large structural sweep. These issues become impor- 
tant in cases IIA-1 and IIA-2 wherein large bending 
deflections occur at the high collective angles, 
and in case IIA-6 where the effects of structural 
sweep are most pronounced. The following 
subsections address these two issues. 
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Torsion Excitation due to Compound Bending 

As given in Reference 1, and as recognized 
elsewhere in the literature, the torsion differen- 
tial equation is comprised of three basic parts. 

The first part consists of the usual elastic 
stiffening terms, and the second consists of 
combinations of distributed moment loadings. The 
third part is comprised of the wholly nonlinear 
torsion loadings accruing from distributed force 
loadings acting on moment arms provided by curvature 
in the elastic axis. As given in Reference 1, the 
torsion equation is given by: 

-/ 

[GJfle + ®>dT + e'eW- EB 2 0e v e''J ^ 


elastic stiffening 

*[- V y " Q *" I » q, d(g) 


moment loadings (15 

+ K'U'k / r ' p » 5 < r 2 )dr 2- / r 'p l5 (^ dr 2 + <W r ' ) ] dr « 

I ' 

- 2 *"'l I,' K f,[ P,^2> dr 2- 

V — -v ' 

curvatures functions of force loadings 


not assured with a small number of modes. This 
accuracy problem is then compounded by the fact that 
the two components of this nonlinear excitation are 
subtractive. This is evidenced by the differencing 
of the section bending stiffnesses as indicated above . 

A second difficulty with using the AEI method 
relates to the assumed space curve character of the 
elastic axis. As such, torsion deflections are 
seen to contribute to inplane and out-of-plane 
deflections in the presence of bending. Thus, an 
analogous nonlinear excitation exists in both the 
flatwise and edgewise bending equations. In the 
framework of the GAOO technology, these nonlinear 
excitations in the bending equations are most 
practically implemented using a ’’force integration" 
approach. Consequently, the use of a AEI mode 
deflection implementation in the torsion equation 
together with a force integration implementation 
in the bending equations results in a (coupled) 
modal mass matrix which is generally nonsymmetric . 

A nonsymmetric mass matrix is not intrinsically 
a weakness for isolated rotor simulation and has 
been successfully used for years in that mode. 
However, the potential exists for spurious diver- 
gent response conditions caused by an inertia 
matrix becoming nonpositive-definite due to this 
deflection dependent nonsymmetry . 

The third difficulty with the Equation 16 for- 
mulation is that it is difficult to include the 
built-in curvature due to structural sweep. 

Equation (16) requires curvature information which 
is not generally available for the built-in 
geometry . 


In Reference 1, the curvatures used in the 
(nonlinear) third portion of the torsion equation 
were assumed to arise entirely from the elastic 
bending deflections, v£ and w^ . As such, it can 
be shown that the nonlinear excitation term in 
Equation (37) can be reduced to a compact expression 
which includes the familiar difference of bending 
stiffness term, AEI ( * EI z -EI y ) : 

{•••} @ = [<u*-Ei y >v;v 

- (e A T + EB 2 (fl' + ^0 e ')S«W'] (16) 

This method for including the effect is 
attractive principally because of its simplicity 
and has been used to good advantage by numerous 
investigators. Three difficulties exist with this 
method of implementation, however. The first 
difficulty relates to the fact that the implementa- 
tion of Equation (16) is based on a "mode deflection" 
description of internal bending moment. The diffi- 
culty with a mode deflection formulation per se is 
two-fold. Studies of the characteristics of 
mode deflection formulations (References 12 and 13) 
have established that convergence to accurate 
representations of internal bending moment is often 


Because of these difficulties, the conven- 
tional AEI approach of Equation (16) was abandoned 
in favor of a "force integration" approach. 
Accordingly, the Galerkin approach is first applied 
to the nonlinear excitation term and then 
integration by parts is used to achieve an inter- 
mediary step needed to eliminate the explicit 
curvature terms: 


+ p y 5 Oo Ye j dr * dr ' + (z * T + V y s dr ' 

My 5 'T-q l5 )/V e z * 


( 17 ) 


Since this term represent the nonlinear 
effects, it is reasonable to use a zeroth order 
approximation to the curvature terms wherein the 
structural sweep is assumed to be "small". With 
this assumption, all the integrals in Equation (17) 
can be evaluated using the deflection correction 
functions defined in Reference 1. 


Thus, Equation (17) becomes: 


•£ r °> {*** l fejlv 05 ® + t v n ®l 

" ^ [ p i 5 cos ® - P^sin©] 

+ f z tf .[ T(w « + z fo EA " Aw (2) - AW ,2) ') + Qy j cos© + q Zj sin©] 

[ T ( y e + yio EA + Av ,2) - Av ( 2 ’)- q^cos® + q y5 sin©]}dr 

(18) 

where : 

r »fij S y 0j(We + z io EA _ Aw-Aw)-(Av ea .-Av ea .) (19a) 

V ^j (v e + yiO EA + Av “ AV) ~ (Aw E Aj + AW ea .) (19b) 

^Ve'+zioE*-^' 2 ’ -AW (2>/ )- (Av ea .- AV ea .') (19c) 

^0j = r 0j ' v e + y,'o eA + Av< 2 >- AV (21 ') - (Aw ea (2), + AW ea ’') 

(19d) 


Equation (18) represents the required form of 
the ’’force integration” implementation of the 
nonlinear torsion excitation term. Upon recog- 
nizing and utilizing various cancellations arising 
in Equation (18) itself and in combination with 
similar terms contained in the © moment loadings 
term, the final most useful form of the torsion 
equation can then be written as: 

fo Ye \ + -J dr= / 0 ‘{ V«5" fy[ p v 5 C0S ® + V in ®l 

+ r *e [ p i, cos ® “ P v sin©] 

] 5 J 5 J 

+ (Aw ea < 2) '+ Aw ea ! 2 >') [t(w, + ... + q^sin©] 

(2)^ (2 ) 1 r 

-(av ea .-av e ;; ) [T(ve+ + Qy 5 cos@Jj dr (20) 


To conclude this subsection, three observations 
can be made of the above formulations: 


1. Equations (19) all reduce to zero for 
zero structural sweep and zero elastic 
deflection, as would be expected from 
the behaviour of Equation (16). 


2. 


In Equation (18), the terms multiplying 
the nonlinear torsion weighting functions 


(r 


yej 


. .) are actually the force and 


moment loadings defined for the linear 
excitations of the bending equations. 

The nonlinear torsion weighting functions, 
Equations (19), thus serve in effect, 
as the virtual deflection functions 
arising from torsion deflections 
appropriate to the bending generalized 
loads. 

3. The validity of the force integration 
approach is enhanced by the fact that 
the resulting terms in the torsion 
equation which represent rows of the 
inertia matrix (reflecting the integra- 
tion of inertia forces) produce 
complete mass matrix symmetry and 
consequently insure positive-definiteness. 

Kinematic Representation for Structural Sweep 

The selected general approach to modeling struc- 
tural sweep is to use the simple well established 
concepts for bending and torsion of straight beams 
as a departure point. Accordingly, blade elastic 
bending is defined by conventional beam bending 
differential equations wherein the usual independent 
spanwise variable is taken to be the arc length 
along the elastic axis. Furthermore, these bending 
differential equations are defined locally using the 
loadings normal to the built-in elastic axis. 

Within this context, explicit elastic bending-torsion 
coupling due to structural sweep is omitted in favor 
of implicit coupling due to inertial, aerodynamic 
and gravitational loadings taken with appropriate 
sweep related kinematics. Within this context, 
the major necessary task in modeling structural 
sweep is to define the kinematics of the blade 
element mass centers and aerodynamic centers as 
explicit functions of the blade modal response 
variables. This subsection addresses this major 
task, from which the formulations of inertial 
aerodynamic and gravity loads follow in a straight- 
forward manner. These subsequent formulations for 
loadings are thus omitted herein for clarity. 

Structural sweep is defined in a general sense 
wherein both inplane and out-of-plane offsets of 
the built-in elastic axis, y 5 ^ and Z 5 , respec- 
tively, are admitted (see Figure 6). The basic 
objectives of the structural sweep related 
reformulations are: ( 1 ) to define a coordinate 

system rotation transformation from the " 5 ” pitch 
axis system to the swept ” 6 ” system (which is 
locally attached to the elastic axis), and ( 2 ) to 
define the deflections in the " 5 ” system as 
functions of the built-in structural sweep and 
the elastic bending and torsion motions, which are 
measured in the " 6 ” system. These two objectives 
must also be met while including the previous G400 
formulations with regard to structural twist. 

The procedure formulated for including these two 
structural elements (sweep and twist) is summarized 
in the material which follows; the reader is 
directed to Reference 4 for a more detailed 
description. 
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The general modeling of the blade and 25 
kinematics due to combined structural twist and 
sweep is accomplished in the following steps: 

1 . The elastic axis of the "equivalent beam" 
described in an above subsection is 
"distorted" back to the original planform 
defined by the built-in structural sweep 
and segment arc length distributions (but 
without pitch or twist). This step essen- 
tially defines the position in space of 
the elastic axis space curve. This posi- 
tioning requires the X 5 , y 5 and Z 5 offset 
distances of the centers of the segments 
as well as projections onto the X5~y5 and 
x 5 -z 5 pl anes of the swept elastic axis 
line segments. These projections define 
the sweep angle distributions, A e ^ and Af,-. 

2 . As shown in Figure 6, the orientations of 
the elastic axis line segments define the 
local "6" coordinate system. x$ is defined 
parallel to the axis of the elastic axis 
line segment; y^ is defined parallel to 
the x^-y5 plane, (+) in leading edge 
direction; is orthogonal to x 6 and y 6 , 
(+) in the normally positive thrusting 
motion. It should be stressed that the 
result of step 1 is to produce, in addi- 
tion to the inplane and out-of-plane 
offsets (Ay 5 and Az 5 ) of the elastic axis 
from the (reference) X5 pitch axis, a 
radial foreshortening ( x^) due to the 
constancy of the total arc length of the 
elastic axis. This Ax 5 foreshortening 

is given by the negative of u e , as 
developed in the previous section. 


3 . The blade segments of the blade configura- 
tion resulting from steps 1 and 2 are then 
pitched and twisted about their respective 
elastic axis line segments (x^ axis) to 
restore the blade back to its original 
built-in, but elastically undeflected 
position. The pitch and twist angles for 
each segment are defined relative to the 
y 6 axis. 


4 . 


The blade is then elastically deflected in 
torsion (0 =£/Y 0 .30, ) about tbe built-in space 

J J J 

curve elastic axis as defined by YlO^ an< ^ 


ziO£A t0 define a first set of "small" incre- 
mental y^ and 25 deflections. 


5 . The blade is then elastically deflected in 
flatwise and edgewise bending (w and v, 
respectively in the presence of the torsion 
deflection) to define a second set of small 
incremental deflections. This second set of 
incremental deflections is measured in the 
"6" coordinate system and is governed by the 
basic G 400 deflection correction transforma- 
tions defined in Reference 1 . 






Fig. 6 Schematics of the " 5 " and "6" coordinate 
systems . 


6. The second set of small incremental "6" 
coordinate system deflections defined in 
step 5 is transformed to the " 5 M coordinate 
system using an Euler angle transformation 
derived from sweep angle projections, 


A e ^ and 


'£5- 


as discussed in above step 1. 


7 . The results of steps 1 , 4 and 6 are then 
combined to define the total ys and Z5 
position vector components. These procedures 
are mathematically described by the following 
material . 


First, the sweep angle projection distributions 
are defined using the built-in elastic axis line 
segment changes per segment length, the (invariant) 
segment arc lengths Ar, together with changes to the 
projection angles caused by elastic torsion 
deflection: 


A V *>" f 15 “ - 4V„; !,I }C0S® 

* UwcY* £ w E» ] 2 ''lS'r^]q 8 J (21) 

A ffi = sin -1 ) cos© 

+ (ZW EA < 2,/ - AV EA ( . 2)/ )sin®] q*. } ( 22) 
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where y 5 ^ and z 5£A are the built-in elastic 
axis offset changes per segment length. For 
consistency with the definitions used for other 
previously defined radial distributions, these 
spanwise variable quantities are considered to be 
’'derived" quantities calculated from the corres- 
ponding quantities defined in the chordwise and 
thicknesswise directions, y 1Q and z 10 , 

respectively. In practice, however, the ^ 

5 coordin 4 te system quantities are the more 
accurately known and the ”10" coordinate system 
quantities are derived using trigonometric reso- 
lution with the local built-in pitch angle. 

The coordinate system transformation relating 
the pitch axis ("5") coordinate system with the 
swept (”6") coordinate system makes use of the 
sweep angle projections given in Equations (21) 
and (22): 


AX fi 


■ M 


Ax, 


(23) 


AX S 


where : 
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[tas h ] ' AXg 

[tas t ] 

AXc ■ 

X 

-sin A e5 

sinA f 5 

sinAe s 

X 

o 

cos Af 5 

cosA f 

T 5 


- X sinAf 5 
cosAf 5 

sinAf 5 sinA, 
cosAf 5 

*5 

— cos Af 5 


(24) 


(25) 
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(28) 

and where v e , w e , Av, Aw, AV, AW are linear and 
nonlinear combinations of q w ., q and q e , as per 
the original G400 structural twist formulations (1). 
Thus, the objectives defined above have been met; 
the addition of structural sweep is accomplished 
while retaining the structural twist formulation. 

The formulation given by Equation (27) together with 
that for radial foreshortening. Equation (11), 
extends the kinematic modeling to applications with 
large structural sweep and moderate structural twist. 
Note that these formulations are generally quite 
nonlinear in the elastic modal response variables, 

^ q 0j - 

Eigensolutions vs. Time-History Solutions 

As shown in Figure 1, the basic G400 mathema- 
tical capability includes both an eigensolution and 
a time-history solution. Yet, despite the known 
advantages of eigensolutions, the time-history 
solution capability was used exclusively and 
produced results which were probably unattainable 
using the conventional eigensolution approach. 

The generally well-identified disadvantages of time- 
history solutions relative to eigensolutions (for 
stability calculations) are: (1) The calculation 

(CPU) time, and hence cost, is at least one order 
of magnitude greater; (2) the calculations 
inherently include the integral order forced 
responses which obscure assessment of the 
transients, and (3) postprocessing is required to 
obtain conventional stability descriptors. 


M- 


0 I 0 
0 0 I 


and where: 

X=%4^sin 2 Ae 5 - sm 2 Af 5 


(26) 


The above development can then be combined to 
yield the required expressions for inplane and 
out-of-plane displacement: 


f y 5 1 f y »EA C0S ®B- Z IO EA 

l z 5 J L y 'O EA sm0B+Zio EA cose 8 J 
+ f <Av EAj" AV ea .)cos ® + (Aw ea . + AW ea .) sin© 

i'l l (Av EAj- ^EAj>sin©-(Aw EA . + AW ea .)cos 0J q ®i 


+ [e][ t as-'] 


) 0 

\ (v e + Av- AV)cos®-(w e - Aw- AW)Sin© > 

V (v e + Av - AVlsin© + (w e - Aw - AW) cos © J 

(27) 


The time-history solution, as formulated and 
implemented in the G400 technology, does not solve 
essentially linearized equations using an appro- 
priate quadrature algorithm. Rather, the dynamic 
equations are retained in their nonlinear 
(implicit) form without recourse to the explicit 
expansion of loadings (as is typically required for 
eigensolutions). For the present study, this 
compact implementation presented clear advantages 
which outweighed the above identified disadvantages: 
(1) the accuracy of the basic physical modeling is 
separated from the issue of selected linearization 
scheme (mathematical modeling); (2) there is no need 
to calculate accurate equilibrium trimmed responses 
(as required for eigensolution linearization 
schemes), and, most significantly, (3) the compact 
implicit modeling scheme allows physical modeling 
modifications to be made easily to the coding and 
then quickly evaluated. It should be stressed that 
these advantages are related mostly to research and 
methodology development issues. For routine produc- 
tion calculations, the cost-effectiveness of eigen- 
solutions is not to be denied. Thus, a synergistic 
relationship is implied between time-history 
solution and eigensolution development. The former 
is the superior physics modeling tool needed by the 
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latter before the mathematical modeling processes of 
linearization should occur. Clearly, the former 
provides an excellent alternate basis for validating 
the latter, whereas, the latter, once validated 
provides superior computational resources to the 

analyst . 

Concluding Remarks 

The challenge posed by the ITR/FRR Methodology 
Assessment study to correlate analyses with detailed 
experimental stability data has borne fruitful 
advances in the development of aeroelastic methodo- 
logy. The United Technologies G400 analysis after 
being upgraded as a result of this study now 
appears to be well validated. Whereas, some of the 
reformulations constituting this upgrading are 
indigenous only to the G400 technology base, others 
appears to have general applicability to the field 
of rotor aeroelastics . These reformulations consti- 
tute, in part, some new solution techniques for some 
old problems: the inclusion of vortex-ring state 

effects into air mass dynamics, the kinematics of 
radial foreshortening, and a more accurate modeling 
of tension. More significantly, these reformula- 
tions constitute solution techniques for the 
relatively new problem area posed by combined 
variable structural sweep and structural twist. 

These latter reformulations should find useful 
application to a wide range of advanced rotor craft, 
such as aeroelastically conformable helicopter rotor 
blades, advanced technology propellers and prop- 
fans. 

References 

1. Bielawa, R. L., "Aeroelastic Analysis for 
Helicopter Rotor Blades with Time-Variable, 
Nonlinear Structural Twist and Multiple 
Structural Redundancy-Mathematical Derivation 
and Program User’s Manual," NASA CR-2638, 
February 1976. 

2. Bielawa, R. L., "A Second Order Nonlinear 
Theory of the Aeroelastic Properties of 
Helicopter Rotor Blades in Forward Flight," 

Ph.D. Thesis, Massachusetts Institute of 
Technology, June 1965. 


3. Bielawa, R. L. , "Aeroelastic Analysis for 
Helicopter Rotors with Blade Appended 
Vibration Absorbers - Mathematical Derivations 
and Program User’s Manual," NASA CR-165896, 

June 1982. 

4. Bielawa, R. L. , Chi, R. M. , Johnson, S. A., and 
Gangwani, S. T. , "Aeroelastic Analysis for 
Propellers - Mathematical Formulations and 
Program User’s Manual," NASA CR-3729, June 1983. 

5. Bielawa, R. L. , Cheney, M. C. Jr., Novak, R. C., 
"Investigation of a Bearingless Helicopter Rotor 
Concept Having a Composite Primary Structure," 
NASA CR-2637 , October 1976. 

6. Sopher, R. , and Kottapalli, S.B.R., "Correlation 
of Predicted Vibrations and Test Data for a Wind 
Tunnel Helicopter Model," 38th Annual National 
Forum Proceedings, May 1982. 

7. Anonymous, "Integrated Technology Rotor 
Methodology Assessment - Request for Proposal, 
NASA Ames Research Center, FRP 2-30221 (BR) , 
September 1980. 

8. Ormiston, R. A., and Bousman, W. G., "A Study of 
Stall-Induced Flap-Lag Instability of Hingeless 
Rotors," American Helicopter Society 29th Annual 
National Forum Proceedings, May 1973. 

9. Peters, D. A., and Gaonkar, G. H. , "Theoretical 
Flap-Lag Damping with Various Dynamic Inflow 
Models," Journal of the American Helicopter 
Society, Vol. 25, No. 3., July 1980. 

10. Payne, P. R. , Helicopter Dynamics and Aero- 
dynamics, Pitman & Sons, Ltd., London, 1959. 

11. Gessow, A., and Myers, G., Aeordynamics of the 
Helicopter, MacMillan Co., 1952, Republished by 
Frederick Ungar Publishing Company, New York, 
1967. 

12. Bisplinghof f , R. H. , Ashley, H. , and Halfman, 

R. L. , Aeroelasticity , Addison-Wesley Publishing 
Co., Inc. 1957. 

13. Bielawa, R. L. , "Blade Stress Calculations Mode 
Deflection vs. Force Integration," Journal of 
the American Helicopter Society, Vol. 24, No. 3, 
July 1978. 


228 



